#include <iostream>
#include <cuda.h>
#include <stdio.h>
#include <IL/il.h>

using namespace std;

#define BLOCK_WIDTH  32

__global__ void kernel_convert_to_grayScale(unsigned char* data, unsigned char* outGreyScale, unsigned int width, unsigned int height){

  	int x = blockIdx.x * blockDim.x + threadIdx.x;
	int y = blockIdx.y * blockDim.y + threadIdx.y;

	if ((y < height) && (x < width)) { 
        
        int idx = (x + width * y) * 3;        // index du pixel courant
 
        outGreyScale[idx/3] = ((data[idx] * 307 + data[idx + 1] * 604 + data[idx + 2] * 113) >> 10);
    }
}

__global__ void kernel_sobel(unsigned char* data, unsigned char* outSobel, unsigned int width, unsigned int height){

  	int x = blockIdx.x * blockDim.x + threadIdx.x;
	int y = blockIdx.y * blockDim.y + threadIdx.y;

	//Sobel Horizontal Mask
	double filtreGx [3][3] = { 

		{1.0,0.0,-1.0},
		{2.0,0.0,-2.0},
		{1.0,0.0,-1.0}
	};
	
		//Sobel Vertical Mask
	double filtreGy [3][3] = { 

		{1.0,2.0,1.0},
		{0.0,0.0,0.0},
		{-1.0,-2.0,-1.0}
	};
	
	double sum = 0.0;
	double valx = 0.0;
	double valy = 0.0;
				
	for(int i = -1; i <= 1; i++) {
				
		for(int j = -1; j <= 1; j++){
			
			int pixel = y*width+x+i*width+j;        // index du pixel courant
			
			if ((y+i) >= 0 && (y+i) < height && (x+j) >= 0 && (x+j) < width){
				
				valx += (double)data[pixel] * filtreGx[ i + 1][ j + 1];
				valy += (double)data[pixel] * filtreGy[ i + 1][ j + 1];
				
			}
			
			//------ recalcule des pixels ----
			
			sum = sqrt(valx*valx + valy*valy);
			
			outSobel[y*width+x] = sum;
		}
	}
}

__global__ void kernel_histo(unsigned char* data, int* outHisto, unsigned int width, unsigned int height){

	int x = blockIdx.x * blockDim.x + threadIdx.x;
	int y = blockIdx.y * blockDim.y + threadIdx.y;

	int idx = (x + width * y);        // index du pixel courant
	
	atomicAdd(&outHisto[data[idx]],1);

}

__global__ void kernel_histoEq(unsigned char* data, unsigned char* outHistoEq, double* histoCumule, unsigned int width, unsigned int height){

	int x = blockIdx.x * blockDim.x + threadIdx.x;
	int y = blockIdx.y * blockDim.y + threadIdx.y;

	int idx = (x + width * y);        // index du pixel courant
	
	outHistoEq[idx] = (256.0  / (double)(width*height)) * histoCumule[data[idx]]-1;

}

/*void initialiseGPU(unsigned char* data, unsigned char* out, int size, int width, int height) {

	cudaEvent_t start, stop0, stop1;
  	float time;

  	cudaEventCreate(&start);
  	cudaEventCreate(&stop0);
  	cudaEventCreate(&stop1);

	unsigned char* dataGPU;
	unsigned char* outGPU;

	cudaEventRecord(start, 0);

	//on alloue sur la mémoire du GPU
	cudaMalloc(&dataGPU, size * sizeof(unsigned char));
	cudaMalloc(&outGPU, size * sizeof(unsigned char));

	// on écrit des 0 dans la zone mémoire alloué pour etre sur qu'il ne récupere pas les données de la derniere exécution 
	cudaMemset(dataGPU, 0, size * sizeof(unsigned char));
	cudaMemset(outGPU, 0, size * sizeof(unsigned char));

	//On recopie dans la mémoire du GPU l'image de départ
	cudaMemcpy(dataGPU, data, size * sizeof(unsigned char), cudaMemcpyHostToDevice);

	//cudaError hr = cudaSuccess; 

	int x = static_cast<int>(ceilf(static_cast<float>(width) / BLOCK_WIDTH));
    int y = static_cast<int>(ceilf(static_cast<float>(height) / BLOCK_WIDTH));
 
    const dim3 grid (x, y);                                // nombre de blocks
    const dim3 block(BLOCK_WIDTH, BLOCK_WIDTH);            // block width: nombre de threads par blocks

	//const dim3 block(24, 24, 1);
 	//const dim3 grid((width/64), (height/64) , 1);

    convert_to_grayScale<<<grid,block>>>(dataGPU, outGPU, width, height);
 	//sobel<<<grid,block>>>(dataGPU, outGPU, width, height);

	cudaEventRecord(stop0, 0);

	cudaEventSynchronize(stop0);

	//hr = cudaDeviceSynchronize();
	//CHECK_CUDA_ERROR(hr, func, "ConvertToGrayScale failed."); 

	// on recopie le résultat dans la mémoire du CPU
	cudaMemcpy(out, outGPU, (width * height) * sizeof(unsigned char), cudaMemcpyDeviceToHost);

	cudaEventRecord(stop1, 0);

	cudaEventSynchronize(stop1);

	cudaEventElapsedTime(&time, start, stop0);
  	cout << "time kernel=" << time << endl;

  	cudaEventElapsedTime(&time, stop0, stop1);
  	cout << "time copy=" << time << endl;

  	cudaEventDestroy(start);
  	cudaEventDestroy(stop0);
  	cudaEventDestroy(stop1);

	cudaFree(dataGPU);
	cudaFree(outGPU); 

}*/

void ConvertToGrayScale(unsigned char* data, unsigned char* outGreyScale, int width, int height) 
{  

	cudaEvent_t start, stop0, stop1;
  	float time;

  	cudaEventCreate(&start);
  	cudaEventCreate(&stop0);
  	cudaEventCreate(&stop1);

	unsigned char* dataGPU;
	unsigned char* outGreyScaleGPU;

	cudaEventRecord(start, 0);


	//on alloue sur la mémoire du GPU
	cudaMalloc(&dataGPU, (width*height*3) * sizeof(unsigned char));
	cudaMalloc(&outGreyScaleGPU, width*height * sizeof(unsigned char));

	cudaMallocHost(&data, (width*height*3) * sizeof(unsigned char));

	// on écrit des 0 dans la zone mémoire alloué pour etre sur qu'il ne récupere pas les données de la derniere exécution 
	cudaMemset(dataGPU, 0, (width*height*3) * sizeof(unsigned char));
	cudaMemset(outGreyScaleGPU, 0, width*height * sizeof(unsigned char));

	for (int i=0; i<5; i++) {

		//On recopie dans la mémoire du GPU l'image de départ
		cudaMemcpyAsync(dataGPU, data, (width * height * 3) * sizeof(unsigned char), cudaMemcpyHostToDevice);

		int x = static_cast<int>(ceilf(static_cast<float>(width/4) / BLOCK_WIDTH));
	    int y = static_cast<int>(ceilf(static_cast<float>(height) / BLOCK_WIDTH));
	 
	    const dim3 grid (x, y);                                // nombre de blocks
	    const dim3 block(BLOCK_WIDTH, BLOCK_WIDTH);            // block width: nombre de threads par blocks

		kernel_convert_to_grayScale<<<grid,block,i>>>(dataGPU, outGreyScaleGPU, width, height);

		//cudaEventRecord(stop0, 0);

		//cudaEventSynchronize(stop0);

		// on recopie le résultat dans la mémoire du CPU
		cudaMemcpyAsync(outGreyScale, outGreyScaleGPU, (width * height) * sizeof(unsigned char), cudaMemcpyDeviceToHost);

	}

	/*cudaEventRecord(stop1, 0);

	cudaEventSynchronize(stop1);

	cudaEventElapsedTime(&time, start, stop0);
  	cout << "time kernel=" << time << endl;

  	cudaEventElapsedTime(&time, stop0, stop1);
  	cout << "time copy=" << time << endl;*/

  	cudaEventDestroy(start);
  	cudaEventDestroy(stop0);
  	cudaEventDestroy(stop1);

	cudaFree(dataGPU);
	cudaFree(outGreyScaleGPU); 

}

void Sobel(unsigned char* data, unsigned char* outSobel, int width, int height)
{

	cudaEvent_t start, stop0, stop1;
  	float time;

  	cudaEventCreate(&start);
  	cudaEventCreate(&stop0);
  	cudaEventCreate(&stop1);

	unsigned char* dataGPU;
	unsigned char* outSobelGPU;

	cudaEventRecord(start, 0);

	//on alloue sur la mémoire du GPU
	cudaMalloc(&dataGPU, (width*height) * sizeof(unsigned char));
	cudaMalloc(&outSobelGPU, width*height * sizeof(unsigned char));

	// on écrit des 0 dans la zone mémoire alloué pour etre sur qu'il ne récupere pas les données de la derniere exécution 
	cudaMemset(dataGPU, 0, (width*height) * sizeof(unsigned char));
	cudaMemset(outSobelGPU, 0, width*height * sizeof(unsigned char));

	//On recopie dans la mémoire du GPU l'image de départ
	cudaMemcpy(dataGPU, data, (width * height) * sizeof(unsigned char), cudaMemcpyHostToDevice);

	int x = static_cast<int>(ceilf(static_cast<float>(width) / BLOCK_WIDTH));
    int y = static_cast<int>(ceilf(static_cast<float>(height) / BLOCK_WIDTH));
 
    const dim3 grid (x, y);                                // nombre de blocks
    const dim3 block(BLOCK_WIDTH, BLOCK_WIDTH);            // block width: nombre de threads par blocks

	kernel_sobel<<<grid,block>>>(dataGPU, outSobelGPU, width, height);

	cudaEventRecord(stop0, 0);

	cudaEventSynchronize(stop0);

	// on recopie le résultat dans la mémoire du CPU
	cudaMemcpy(outSobel, outSobelGPU, (width * height) * sizeof(unsigned char), cudaMemcpyDeviceToHost);

	cudaEventRecord(stop1, 0);

	cudaEventSynchronize(stop1);

	cudaEventElapsedTime(&time, start, stop0);
  	cout << "time kernel=" << time << endl;

  	cudaEventElapsedTime(&time, stop0, stop1);
  	cout << "time copy=" << time << endl;

  	cudaEventDestroy(start);
  	cudaEventDestroy(stop0);
  	cudaEventDestroy(stop1);

	cudaFree(dataGPU);
	cudaFree(outSobelGPU); 
}

void Histo(unsigned char* data, int* histo, int width, int height)
{

	cudaEvent_t start, stop0, stop1;
  	float time;

  	cudaEventCreate(&start);
  	cudaEventCreate(&stop0);
  	cudaEventCreate(&stop1);

	unsigned char* dataGPU;
	int* outHistoGPU;

	cudaEventRecord(start, 0);

	//on alloue sur la mémoire du GPU
	cudaMalloc(&dataGPU, (width*height) * sizeof(unsigned char));
	cudaMalloc(&outHistoGPU, 256 * sizeof(int));

	//On recopie dans la mémoire du GPU l'image de départ
	cudaMemcpy(dataGPU, data, (width * height) * sizeof(unsigned char), cudaMemcpyHostToDevice);
	cudaMemcpy(outHistoGPU, histo, 256 * sizeof(int), cudaMemcpyHostToDevice);	

	int x = static_cast<int>(ceilf(static_cast<float>(width) / BLOCK_WIDTH));
    int y = static_cast<int>(ceilf(static_cast<float>(height) / BLOCK_WIDTH));
 
    const dim3 grid (x, y);                                // nombre de blocks
    const dim3 block(BLOCK_WIDTH, BLOCK_WIDTH);            // block width: nombre de threads par blocks

	kernel_histo<<<grid,block>>>(dataGPU, outHistoGPU, width, height);

	cudaEventRecord(stop0, 0);

	cudaEventSynchronize(stop0);

	// on recopie le résultat dans la mémoire du CPU
	cudaMemcpy(histo, outHistoGPU, 256 * sizeof(int), cudaMemcpyDeviceToHost);

	cudaEventRecord(stop1, 0);

	cudaEventSynchronize(stop1);

	cudaEventElapsedTime(&time, start, stop0);
  	cout << "time kernel=" << time << endl;

  	cudaEventElapsedTime(&time, stop0, stop1);
  	cout << "time copy=" << time << endl;

  	cudaEventDestroy(start);
  	cudaEventDestroy(stop0);
  	cudaEventDestroy(stop1);

	cudaFree(dataGPU);
	cudaFree(outHistoGPU); 
}

void HistoEq(unsigned char* data, unsigned char* outHistoEq, double* histoCumule, int width, int height) {

	cudaEvent_t start, stop0, stop1;
  	float time;

  	cudaEventCreate(&start);
  	cudaEventCreate(&stop0);
  	cudaEventCreate(&stop1);

	unsigned char* dataGPU;
	unsigned char* outHistoEqGPU;
	double* histoCumuleGPU;

	cudaEventRecord(start, 0);

	//on alloue sur la mémoire du GPU
	cudaMalloc(&dataGPU, (width*height) * sizeof(unsigned char));
	cudaMalloc(&outHistoEqGPU, (width*height) * sizeof(unsigned char));
	cudaMalloc(&histoCumuleGPU, 256 * sizeof(double));

	// on écrit des 0 dans la zone mémoire alloué pour etre sur qu'il ne récupere pas les données de la derniere exécution 
	cudaMemset(dataGPU, 0, (width*height) * sizeof(unsigned char));
	cudaMemset(outHistoEqGPU, 0, (width*height) * sizeof(unsigned char));

	//On recopie dans la mémoire du GPU l'image de départ
	cudaMemcpy(dataGPU, data, (width * height) * sizeof(unsigned char), cudaMemcpyHostToDevice);
	cudaMemcpy(histoCumuleGPU, histoCumule, 256 * sizeof(double), cudaMemcpyHostToDevice);

	int x = static_cast<int>(ceilf(static_cast<float>(width) / BLOCK_WIDTH));
    int y = static_cast<int>(ceilf(static_cast<float>(height) / BLOCK_WIDTH));
 
    const dim3 grid (x, y);                                // nombre de blocks
    const dim3 block(BLOCK_WIDTH, BLOCK_WIDTH);            // block width: nombre de threads par blocks

	kernel_histoEq<<<grid,block>>>(dataGPU, outHistoEqGPU, histoCumuleGPU, width, height);

	cudaEventRecord(stop0, 0);

	cudaEventSynchronize(stop0);

	// on recopie le résultat dans la mémoire du CPU
	cudaMemcpy(outHistoEq, outHistoEqGPU, (width * height) * sizeof(unsigned char), cudaMemcpyDeviceToHost);

	cudaEventRecord(stop1, 0);

	cudaEventSynchronize(stop1);

	cudaEventElapsedTime(&time, start, stop0);
  	cout << "time kernel=" << time << endl;

  	cudaEventElapsedTime(&time, stop0, stop1);
  	cout << "time copy=" << time << endl;

  	cudaEventDestroy(start);
  	cudaEventDestroy(stop0);
  	cudaEventDestroy(stop1);

	cudaFree(dataGPU);
	cudaFree(outHistoEqGPU); 
	cudaFree(histoCumuleGPU);

}


int main(int argc, char **argv) {


	unsigned int image;

	// Initialisation de la bibliothèque
	ilInit();
	ilEnable(IL_ORIGIN_SET);

	// Création d'une image
	ilGenImages(1, &image);

	// Activation de l'image
	ilBindImage(image);

	// Placement des données du fichier passer en paramètre, dans l'image courante
	if (argc > 1) {
		ilLoadImage(argv[1]);
	}
	else 
    {
        cerr << "Veuillez donner une image à traiter." << endl; 
        exit(1); 
    } 

	int width, height, bpp, format;

	// Récupération des informations sur l'image : largeur, hauteur, octets par pixel, format
	width = ilGetInteger(IL_IMAGE_WIDTH);
	height = ilGetInteger(IL_IMAGE_HEIGHT); 
	bpp = ilGetInteger(IL_IMAGE_BYTES_PER_PIXEL);
	format = ilGetInteger(IL_IMAGE_FORMAT);

	cout<<"largeur : "<<width<<endl;
	cout<<"hauteur : "<<height<<endl;
	cout<<"bpp : "<<bpp<<endl;
	cout<<"format : "<<format<<endl;

	// Récupération du contenu de l'image (pixels)
	unsigned char* data = ilGetData();

	// Création des images de sortie
	unsigned char* outGreyScale = (unsigned char*)malloc(width * height);
	unsigned char* outSobel = (unsigned char*)malloc(width * height);
	unsigned char* outHisto = (unsigned char*)malloc(width * height);

	int* histo = new int[256];
    double* histoCumule = new double[256];
	double somme = 0.;

	// initialisation du tableau de l'histogramme
	for (int i=0; i<256; i++) {
	
		histo[i] = 0;
	}

	try 
    { 
    	//initialiseGPU(data, outGreyScale, width, height);
        ConvertToGrayScale(data, outGreyScale, width, height);
        Sobel(outGreyScale, outSobel, width, height);
        Histo(outGreyScale, histo, width, height);
    }
    catch(exception& e) 
    { 
        cerr << endl << "ERROR: " << e.what() << endl; 
        exit(1); 
    }

    // calcul de l'histogramme cumulé
	for (int i=0; i<256; i++) {
		
		somme += histo[i];
		histoCumule[i] = somme;

	}

	try 
    {
        HistoEq(outGreyScale, outHisto, histoCumule, width, height);
    }
    catch(exception& e) 
    { 
        cerr << endl << "ERROR: " << e.what() << endl; 
        exit(1); 
    }


    //############### ECRITURE IMAGES ###############
	
	//convertion de l'image en valeur luminence
	//ilConvertImage(IL_LUMINANCE, IL_UNSIGNED_BYTE);
	ilTexImage(width, height,0,1,IL_LUMINANCE, IL_UNSIGNED_BYTE, outGreyScale);
	
	// Le contenu de l'image de sortie est dans le tableau out
	ilSetData(outGreyScale);

	// Si le fichier existe déjà , écrase son contenu
	ilEnable(IL_FILE_OVERWRITE);

	// Sauvegarde de l'image
	ilSaveImage("outGreyScale.png");

	ilDeleteImages(1, &image); 

	// Image avec Sobel
	ilTexImage(width, height,0,1,IL_LUMINANCE, IL_UNSIGNED_BYTE,outSobel);
	ilSetData(outSobel);
	ilEnable(IL_FILE_OVERWRITE);
	ilSaveImage("outSobel.png");

	// Image pour l'histo
	ilTexImage(width, height,0,1,IL_LUMINANCE, IL_UNSIGNED_BYTE,outHisto);
	ilSetData(outHisto);
	ilEnable(IL_FILE_OVERWRITE);
	ilSaveImage("outHisto.png");

	free(outGreyScale);
	free(outSobel);
	free(outHisto);
 
    cout << "Done!" << endl << endl; 
    
    return 0;
}


/*

32 block

nbBlock2 = height mod 32

nbBlock1 = 32 - nbBlock2

nbLigneBlock1 height / 32

nbLigneBlock2 = (height / 32) + 1







nbLignesParBlock

nbPixelParBlock = nbLignesParBlock * width

1024 threads par block

nbThreads2 = nbPixelBlock mod 1024

nbThreads1 = 1024 - tmp

nbPixelThread1 = nbPixelParBlock / 1024

nbPixelThread1 = (nbPixelParBlock / 1024) + 1


	cudaDeviceProp prop;
	cudaGetDeviceProperties(&prop,0);
	int w = prop.warpSize;
	int h = prop.maxThreadsPerBlock / w;

	cout<<w<<" "<<h<<endl;
	


nvcc Prog_Cuda_chargement_image.cu -lIL -O3 -arch=sm_12
*/